前回は線形計画法(Linear Programming)について学びました。
整数計画法(Integer Programming)はそれによく似ているのですがもっと応用範囲が広く、
解くのが難しいものです。
整数計画法とは、線形計画法の変数に整数を含むものをいいます。
変数が全て整数変数の全整数計画(Integer Programming)と、
変数の一部が整数変数の混合整数計画(Mixed Integer Programming)があります。
前者は略してIP、後者はMIPとよく略されます。
ここでは特に区別せずに扱うものとします。
さて、皆さん整数と聞くとどういうイメージがあるでしょうか？
小学校で初めて習うのが自然数(正の整数)ですから、簡単だと思うかもしれません。
しかし、それは大きな間違いです。
MIPはLPと比較して数百倍も困難な問題で、下手をすると問題のサイズが大きいと
宇宙が終わるまで計算しても解けません。
数十万の変数でも扱えるLPと比べて圧倒的に難しい問題になります。
GLPKは厳密解法(exact method)といって、最適解を計算しますが、
それでは数百変数もあれば限界だと思ってください。
なお近似解法(heuristics method, approximate method)といって、
厳密解ではないがそれほど悪くない解を求める方法もあります。
それはプログラミングを要するので、また次の機会に任せることにしましょう。
いろいろ書いてきましたが、GLPKではLPとMIPの違いをそれほど考えることはありません。 ただ気をつけなければいけないのは、それほど変数の数を多くできないことです。 さっそくですが、問題例を見ていきましょう。
ある外国での話です。
この国では£1, £10, £25の3種類の硬貨が発行されています。
このとき、£44のおつりを返すにはどうしたら
枚数を最小にできるでしょうか？
これは整数計画の一例です。
直感的には、一番大きな硬貨から割り当てていけばうまくいきそうな気がします。
すると、
£25 ・・・ 1枚 （残り£19)
£10 ・・・ 1枚 (残り£9)
£1 ・・・ 9枚 合計£44, コイン数合計11枚
このような方法を欲張り法(greedy method)といい、近似解法の1種です。
つまり、よいと考えられる方法からとりあえず選んでいくわけです。
しかしちょっと考えれば分かるのですが、
本当の最適解は
£25 ・・・ 0枚 （残り£44)
£10 ・・・ 4枚 (残り£4)
£1 ・・・ 4枚 合計£44, コイン数合計8枚
となります。
この例から分かるように欲張り法は必ずしもうまくいかない
ことが分かると思います。
ちなみに日本の硬貨は、欲張り法の解と、最適解が一致するように発行されていますので
ご安心を。
では、これをGLPKを使って解いてみましょう。
まずは集合を考えるのでしたね。
この問題ではコインが重要なので、集合に置いてみましょう。
この集合をCOINとします。
次にパラメータです。
関係あるのは、「コインの価値」と「合計金額」でしょうか。
それぞれ、value, totalとしましょう。
valueはコインによって違うので、value{COIN}とするのでした。
次は変数です。
目的は枚数の最小化なので、それぞれのコインの枚数を変数としましょう。
これをX{COIN}とします。
目的関数は、コインの枚数の総和を最小化するとすればよいでしょう。 制約条件は、合計金額がおつりと一致すればよいですね。
これを数式で表すと以下のようになります。
これをGLPKで書き直してみましょう。
#coin.mod
set COIN;
param value{COIN};
param total;
var X{COIN} ,>=0 ,integer;
minimize Number: sum{c in COIN}X[c];
s.t. Total: sum{c in COIN}value[c]*X[c] = total;
#coin.dat set COIN:= 1 10 25; param value:= 1 1 10 10 25 25; param total:=44;
注目すべきなのは整数の欄です。
var X{COIN} ,>=0, integer;
となっています。
intergerは新しいキーワードで、ご想像通り、整数変数であることを表します。
あとはglpsol -m coin.mod -d coin.dat -o coin.sol
と解いてやれば、結果が得られます。
Problem: coin
Rows: 2
Columns: 3 (3 integer, 0 binary)
Non-zeros: 6
Status: INTEGER OPTIMAL
Objective: Number = 8 (MINimum) 1.76 (LP)
No. Row name Activity Lower bound Upper bound
------ ------------ ------------- ------------- -------------
1 Number 8
2 Total 44 44 =
No. Column name Activity Lower bound Upper bound
------ ------------ ------------- ------------- -------------
1 X[1] * 4 0
2 X[10] * 4 0
3 X[25] * 0 0
End of output
Objective: Number = 8 (MINimum) 1.76 (LP)
というのが最適値です。
8は今求めたかった整数の最適値です。
もう一つの1.76 (LP)というのは、線形緩和解です。
Columnの部分のActivityを見ると最適解は先ほど説明したように、
£1を4枚、£10を4枚で計8枚であることも分かります。
ここでは詳しくふれませんが、内部では分枝限定法(branch and bound method)という
厳密解法を用いており、その手がかりとして整数条件を無視(緩和)することによって
一度LPとして解いた値を利用します。それが1.76という値です。
これ、どこかで見ませんでしたか？
実は、先ほど補足で44/25=1.76とした近似解と等しいのです。
要はこの場合のLPは割り算で十分なほど簡単だということですね。
IPやMIPでは特に0か1しかとらない整数変数があり、これを特に
0-1変数(0-1 variable, binary variable)といいます。
これを含むIPを0-1整数計画にといいます。
0-1変数はやるかやらないかといった2択の意志決定を表すのに役立ち
整数計画の中でもっとも多く使われます。
IPは計算量が大きいですから、単純に変数が整数だからといって
先ほどのように多値の整数計画はあまり使うべきではありません。
例えば、150.1を150に丸めてもそれほど問題はないかもしれません。
このような場合、LPで解いて適当な値に丸めるべきです。
もちろん、厳密性を求められる場合はIPで解いてください。
コインの例の場合、おつりを間違えたらお客さんは猛烈に怒り狂うでしょうし
絶対値が小さいので、例えば4.1を4に近似したら誤差が大きく、
全然違う結果になりかねません。
状況に応じて臨機応変に対応しましょう。
あるところにクマさんのぬいぐるみ専門のプロどろぼうがいました。
彼が持っている袋は特別防護膜つきの貴重品ですが、最大7kgしか入りません。
クマさんはそれぞれ以下のようになっています。
| 重さ | 価値 | |
| A(極小) | 2kg | 16万ドル | B(小) | 3kg | 19万ドル | C(中) | 4kg | 23万ドル | D(大) | 5kg | 28万ドル |
「数理計画とは？」と同じように
こんなときは、1kg当たりの価値が最大のものからとっていくのが自然だと思います。
すると
A(極小) ・・・ 16/2 = 8万ドル/kg
B(小) ・・・ 19/3 = 6.3万ドル/kg
C(中) ・・・ 23/4 = 5.8万ドル/kg
D(大) ・・・ 28/5 = 5.6万ドル/kg
ですから、価値の高い順に、最初にA(2kg)を、次にB(3kg)をとります。
すると合計5kgになり、これ以上入らないので整数解35万ドルです。最適でしょうか？
これも欲張り法は通用しません。
最適解は、A(2kg)とD(5kg)をとって、合計7kgちょうど、最適解44万ドルです。
欲張り法の失敗原因は、残り2kgの隙間ができてしまったからです。
補足:
ちなみに、実数解では欲張り法で上界(upper bound)が分かります。
この場合実数解とは、やはり価値の高い順からとっていき、
A(2kg)とD(5kg)のあと、残り2kgの隙間には、C(4kg)のクマさんを半分に切って
入れるということです。
よって、16+19+23/2 = 46.5万ドルとなります。
これは、クマさんを傷つけてはいけない(整数条件)を無視(緩和)して
得られた解です。
条件を緩くしているので、整数の最適解より必ず同じか大きいことが保証されています。
これを上界といいます。
2のコインの例では、線形緩和により下界(lower bound)が分かったということです。
下界は最適解より小さいと保証されている値です。
GLPKが分枝限定法でLPを解くのは、線形緩和で良い上界(or下界)が得られるからです。
定式化はこれまでと同じように簡単です。
まず集合は、クマさんが必要なのでBearとします。
パラメータは、まずクマさん(Bear)によって変わるものとして
重さと価値あります。これをそれぞれweight、valueとします。
また、袋の重量サイズが必要なので、sizeとします。
変数はどれを取るか？という0-1変数が必要です。
これをXとしておきます。
クマさんをとったら1、とらなかったら0と考えてください。
目的関数は、とったクマさんの価値の合計です。
とったクマさんの変数は1になり、とらなかったクマさんは0になるので
価値×変数がとったクマさんの価値になります。
なぜなら、とらなかったクマさんは変数が0なので、かけざんした価値は0です。
よって、この総和が目的関数になります。
制約条件は、とったクマさんの重量の合計が、袋のサイズ以下であることです。
目的関数と同様の考え方で、重量×変数の総和がサイズ以下であることが条件です。
これらを数式で書くと以下のようになります。
これをGLPKで書くと以下のようになります。
#knap.mod
set Bear;
param weight{Bear};
param value{Bear};
param size;
var X{Bear} binary;
maximize Profit: sum{b in Bear}value[b]*X[b];
s.t. Weight: sum{b in Bear}weight[b]*X[b] <= size;
#knap.dat set Bear:= a b c d; param: weight value:= a 2 16 b 3 19 c 4 23 d 5 28; param size:=7;
注目なのは、変数のところの
var X{Bear} binary;
binaryは新しいキーワードで、変数が0-1変数であることを
表します。
これをGLPKで解けば、結果が得られます。
Problem: nap
Rows: 2
Columns: 4 (4 integer, 4 binary)
Non-zeros: 8
Status: INTEGER OPTIMAL
Objective: Profit = 44 (MAXimum) 46.5 (LP)
No. Row name Activity Lower bound Upper bound
------ ------------ ------------- ------------- -------------
1 Profit 44
2 Weight 7 7
No. Column name Activity Lower bound Upper bound
------ ------------ ------------- ------------- -------------
1 X[a] * 1 0 1
2 X[b] * 0 0 1
3 X[c] * 0 0 1
4 X[d] * 1 0 1
End of output
結果は、
Objective: Profit = 44 (MAXimum) 46.5 (LP)
とあるので、最適値が44万ドルであり、線形緩和解(上界)が46.5万ドルであることが分かります。
さきほどの結果と一致しますね。
また、最適解はX[a]とX[d]のActivityが1になっているので、AとD、つまり極小と大のクマさんを
盗めばいいのですね。
今回は整数計画法について学んできました。
GLPKにはだいぶ慣れてきたでしょうか？
これを解くのは(計算量的に)難しいのですが
GLPKを使えばLPと同じようにモデリングできることも分かったと思います。
IPやMIPにはまだまだ応用例があって、固定費用付き問題や
論理条件なども考えることができます。
次回は整数計画の応用例をいくつかふれていきたいと思います。
それでは！